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ABSTRACT 

A new method of accurately calculating the capacitance of realistic ice particles is described: such values 
are key to accurate estimates of deposition and evaporation (sublimation) rates in numerical weather models. 
The trajectories of diffusing water molecules are directly sampled, using random 'walkers'. By counting 
how many of these trajectories intersect the surface of the ice particle (which may be any shape) and how 
many escape outside a spherical boundary far from the particle, the capacitance of a number of model ice 
particle habits have been estimated, including hexagonal columns and plates, 'scalene' columns and plates, 
bullets, bullet-rosettes, dendrites, and realistic aggregate snowflakes. For ice particles with sharp edges and 
corners this method is an efficient and straightforward way of solving Laplace's equation for the capacitance. 
Provided that a large enough number of random walkers are used to sample the particle geometry (~ 10^) the 
authors expect the calculated capacitances to be accurate to within ~ 1%. The capacitance for our modelled 
aggregate snowflakes (C'/Dmax = 0.25, normalised by the maximum dimension Dmax) is shown to be in 
close agreement with recent aircraft measurements of snowflake sublimation rates. This result shows that 
the capacitance of a sphere {C/ Dmax = 0.5) which is commonly used in numerical models, overestimates the 
evaporation rate of snowflakes by a factor of two. 

The effect of vapor 'screening' by crystals growing in the vicinity of one another has also been investigated. 
The results clearly show that neighbouring crystals growing on a filament in cloud chamber experiments can 
strongly constrict the vapor supply to each other, and the resulting growth rate measurements may severely 
underestimate the rate for a single crystal in isolation (by a factor of 3 in our model setup). 



1. Introduction 

The growth and evaporation of ice particles by 
diffusion of water vapor on to and away from the 
surface of ice particles are fundamental processes for 
the development of ice clouds and precipitation. The 
density of water vapor p around a (stationary) ice 
particle is governed by Laplace's equation: 



(1) 



under steady-state conditions (Pruppacher and Klett 
1997). Typically the separation between ice particles 
in a cloud is such that each ice particle may be consid- 
ered in isolation against a background vapor density 
Paa ■ If the vapor density at the surface Ps is assumed 
to be constant, this leads to the growth rate: 



D Vp-ds^ AirDCip^ - p,), (2) 



dm 
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by application of Gauss's law over the particle surface 
s, where D is the diffusion coefficient for water vapor 
in air, and the capacitance C characterizes the shape 
and size of the ice particle. 

The assumption of constant p^ implies a constant 
surface temperature. The temperature at any 
given part of the surface is determined by both the 
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flux of vapor on to that part of the crystal, resulting 
in the release of latent heat; and by the rate at which 
that heat is conducted away from the crystal. Since 
the equations for diffusion of vapor and transport of 
heat take the same form, even if a portion of the sur- 
face has a larger flux of vapor on to it and becomes 
hotter through the release of latent heat, it should 
lose that extra heat at a correspondingly higher rate, 
resulting in a constant temperature across the whole 
particle surface. In practice the thermal conductiv- 
ity of ice is approximately two orders of magnitude 
larger than the thermal conductivity of moist air, so 
that even if there is some surface migration of the wa- 
ter molecules, the surface temperature should remain 
almost uniform. 

The assumption of constant p^^T^ motivates an 
analogy with electrostatics, where the results for the 
capacitance C are well known for simple geometries 
(see table 1). However, exact solutions for the non- 
smooth shapes of natural ice particles are not avail- 
able, and in numerical models (eg. Wilson and Bal- 
lard 1999) a capacitance based on one of the shapes 
in table 1 is usually applied in its place. It is far from 
clear whether this approximation is a reasonable one, 
or what the capacitance of realistic ice particles actu- 
ally is, particularly for large particles which may be 
complex aggregates. Chiruta and Wang (2003) sug- 
gested a reflnement to this situation by approximat- 
ing bullet rosette crystals by a set of smooth lobes and 
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Table 1 : Theoretical capacitances for simple shapes (McDon- 
ald 1963) . Note that electrical capacitances (in Farads) axe 
usually normalised as iweC, where e is the permittivity of the 
surrounding medium in Fm~^. 
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Capacitance 



Sphere, radius r 
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= r 






Thin disc, radius r 
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Prolate spheroid: major. 
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Oblate spheroid: major. 
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solving equation 1 numerically using a finite-element 
method to obtain the vapor density p around the par- 
ticle, and capacitance C. Recently, the same authors 
have also applied this procedure to solid and hollow 
columns, estimating the capacitance for five model 
aspect ratios (Chinita and Wang 2005). McDon- 
ald (1963) and Podzimek (1966) used metal models 
to simulate realistic ice particles and measure their 
capacitance; however experimental uncertainties lim- 
ited the accuracy of such estimates. 

In this paper a Monte Carlo method for calcu- 
lating the capacitance of realistic ice particles is de- 
scribed. The trajectories of diffusing water molecules 
onto the surface of the ice particle are directly sam- 
pled using random walks. The fraction of walks which 
intersect the modelled ice particle provide an estimate 
for the fiux of water molecules onto the particle, and 
therefore its capacitance. For the non-smooth shapes 
of realistic ice particles, this sampling approach turns 
out to be an accurate and efficient method of solving 
Laplace's equation for C. Also, since we are sam- 
pling steady-state diffusion onto a stationary crystal, 
statistics for C can be built up sampling one random 
walk trajectory at a time, and this means that very 
little computer memory is required. It is interesting 
to note that random walker sampling has recently 
been applied in the electrostatics community to calcu- 
late electrical capacitances for conductors with sharp 
edges (eg. Hwang and Mascagni 2004) since it by- 
passes many of the artefacts introduced by boundary- 
element, finite-difference and finite-element methods 
which can cause systematic errors (Wintle 2004). 

The capacitance of a number of model pristine 
ice crystal types is sampled using this method, along 
with the capacitance of some aggregate snowflakcs 
from the simulations of Westbrook et al (2004). 
These shapes are much more realistic than those that 
are used in numerical models at present, and there- 
fore we should have much more confidence that the 
results given here accurately represent the capaci- 
tance of natural ice particles. This is very important 



if quantitative comparisons are to be made between 
observations and model predictions of deposition and 
evaporation (eg. Forbes and Hogan 2006). 

2. Method 

The simplest Monte Carlo model for solving (1) 
runs as follows. A random walker is placed on a 
large sphere, radius R^a- Its position is chosen at 
random with uniform probability over the sphere's 
surface. This represents the far-field where the sur- 
faces of equal vapor density are spherical. Once re- 
leased from this sphere, the walker is allowed to dif- 
fuse around, taking steps much smaller than any char- 
acteristic length scales in the target ice particle (say 
0.1% of the smallest length scale for accurate sam- 
pling). Each step is taken in a random direction. 
The walker's motion is tracked until either i) it hits 
the surface of the particle, or ii) it escapes outside 
the large sphere (where the 'background' vapor 
density p^^ is fixed). Note that we have treated the 
ice particle as a perfect absorber with p^ = 0, but 
from (2) we see that our calculated capacitances will 
be valid for any value of p^ (including both depo- 
sition and evaporation) provided that it is constant 
over the ice particle surface. The process is repeated 
for a large number of walkers until an accurate sam- 
ple of hits and misses has been obtained: then the 
capacitance may be calculated as C = / x R^ where 
/ is the fraction of walkers that hit the ice parti- 
cle. It is emphasized that the values of the capaci- 
tance C estimated using the above method are the 
same irrespective of the particle surface vapor den- 
sity/temperature: the capacitance is a function of the 
shape and size of the particle alone, provided that p^ 
and Ts are constant. 

In practice, the method outlined above is very in- 
efficient, with a large amount of time being spent 
tracking the walkers as they take very small steps 
through large regions of empty space. By making 
use of the isotropic nature of random walks however, 
a much more efficient algorithm can be constructed. 
The first improvement is to start the walkers uni- 
formly on a much closer sphere which just encloses 
the ice particle (radius R - see figure 1). The distri- 
bution of walkers arriving from R^ onto this closer 
spherical surface for the first time is uniform (because 
there are no sources or sinks of vapor outside R to 
distort the spherical symmetry), and therefore start- 
ing the walks on this much closer sphere introduces 
no bias at all on the statistics. The only caveat is 
that the walkers must still be tracked all the way out 
to i?oo before they can be assumed to have escaped 
to infinity. It is important to emphasise that starting 
the walkers uniformly on R does not force the total 
vapor density to be uniform on R (which would be 
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Fig. 1: Diagram illustrating the random walker method. 
Walkers start from the sphere R enclosing the particle, and 
are tracked along a random walk until they i) hit the particle, 
or ii) escape beyond a distant spherical outer boundary i?oo- 
Note that outside R the walker is able to take big jumps with- 
out any change in its statistics; inside R, the walkers can take 
a jump equal to the distance to the closest point on the sur- 
face of the ice particle. In both cases the jumps are in random 
directions. 



incorrect). We simply make use of the statistical fact 
that walkers arriving from R^^ and passing through 
the surface R for the first time do so randomly and 
uniformly. The total vapor density on R is very much 
non-uniform, and is distorted by the non-spherical ice 
particle, which asserts its shape through the statistics 
of which of the trajectories that pass through R ac- 
tually hit the particle directly, and which trajectories 
escape to pass through the surface of R again (and 
perhaps several times) before hitting the particle or 
escaping to infinity. 

As discussed above, although each walker starts 
on R it must be tracked back to R^ before it is 
assumed to have escaped to infinity (otherwise we 
would be forcing the vapor density to be uniform on 
R rather than on i?oo). Tracking walkers as they move 
outside of R is not a great burden however: the walk- 
ers may take large steps (as large as the shortest dis- 
tance back to the sphere R) in a random direction 
since there are no sources or sinks of vapor within 
that distance and therefore the probability distribu- 
tion for the new position of the walker is isotropic 
(and would be even if the walker took a series of very 
small steps in random directions to travel the same 
linear distance). 

While the random walker is inside R we must be 
more careful. We can still optimize the step length 



however, by calculating the minimum distance from 
the walker's current position to the closest point on 
the ice particle, and setting this distance to be the 
length of the next step. The walker is assumed to 
have hit the ice particle if it comes within some small 
distance 6 of one of the surfaces. The improvements 
described above result in a much more efficient al- 
gorithm, and importantly no extra approximations 
have been made. As before, the process is repeated 
for a large number of walkers so that the particle's 
geometry is accurately sampled, and the capacitance 
is given hy C = f x R. 

The accuracy of the method is determined by 
three factors: i) the number of walkers used; ii) how 
far away the 'infinite' boundary R^ is placed; and 
iii) the thickness of the thin absorbing layer S. In 
the calculations presented here the outer boundary 
is placed at R^o = 500i?, and the absorption layer 
thickness S is set to be less than 0.1% of the smallest 
side length of any of the faces on the particle. There 
are alternative methods that bypass the construction 
of such a layer (Mascagni and Hwang 2003); how- 
ever, for our purposes, simply using a small value of 
6 should be sufficient. To test the accuracy of the 
method, the capacitance of a unit cube has been cal- 
culated, which is known to within a tolerance of 10"'' 
(Hwang and Mascagni 2004). Provided a large num- 
ber of walkers are used (> 10**), the calculation using 
the method described above shows excellent agree- 
ment to well within 1% of this value, as shown in fig- 
ure 2. Sensitivity tests indicate that increasing R^ 
to a value of lOOOi? has no effect on the estimated 
capacitance to within ~ 0.1%*. 

The algorithm has also been tested against par- 
ticles where the analytic solution is known (see ta- 
ble 1). For a thin circular disc with unit radius and 
thickness of 0.001, the capacitance is estimated to be 
C = 0.639; which matches the exact analytical re- 
sult for an infinitely thin disc (C = 2/n = 0.637) to 
within 1%. For a unit sphere we find that our method 
also agrees with the theoretical value (C = 1) to well 
within 1%. In both cases 10" walkers were used. 



3. Results 

Having established that the method is accurate, 
we proceed to apply it to a number of model ice par- 
ticle habits, and the results of this are described be- 
low. 



*We note that it is in fact possible to use the Green's func- 
tion for a point charge outside a grounded sphere to effectively 
place the outer boundary at infinity, removing this source of 
error entirely. This approach also removes the need to track 
the walkers outside of R: see Zhou et al (1994) for details. 
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Fig. 2: Capacitance of unit cube as sampled by an increasing 
number of random walkers. Dashed line shows the theoreti- 
cal value of 0.6607 (to 4 d.p.) as calculated by Hwang and 
Mascagni (2004). The outer boundary was set at Roo = 500R, 
and the thickness of the absorption layer was <5 = 0.001 (0.1% 
of the side length). 

a. Hexagonal columns and plates 

Close to cloud top, one of the most common pris- 
tine ice crystal habits is the hexagonal prism (Prup- 
pacher and Klett 1997). Here the random walker 
method has been used to calculate the capacitance of 
these shapes, and two examples of this kind of crys- 
tal (a column and a plate) are shown in figure 3. We 
define the width 2o to be the maximum span across 
the basal (hexagonal) crystal face, and the length L 
as the span of the crystal perpendicular to the basal 
face. The aspect ratio of the crystal is defined as 
A = L/2a, ie. columns correspond to ^ > 1, plates 
to ^ < 1. 

The capacitance of each crystal was sampled 
by 250,000 walkers, each simulation taking approx- 
imately three minutes on a typical desktop PC. The 
convergence of the sampled capacitances with the 
number of walkers used was consistent with that for 
a cube, and as a result we estimate that the results 
given here for C are accurate to within 1%. The 
distribution of random walkers incident on the faces 
of a column are shown in side projection in figure 4: 
as one would expect from equation 1, walkers impact 
on the particle all across the surface, but are strongly 
concentrated around the particle edges which show 
up as strong dark lines on the figure, reflecting the 
high flux of vapor | Vp| on to these regions. 

Figure 5 shows the capacitance of hexagonal 
prisms C'c ilculated for a constant width (a = 1), where 
the aspect ratio is varied between 0.01 and 50. The 
equation: 

C = 0.58(l+0.95^"-"5)a (3) 




Fig. 3: Geometry of pristine hexagonal column and plate ice 
particles. Width 2a defined as maximum span across basal 
(hexagonal) face; length L is the thickness of the plate, or 
length of the column (ie. the maximum span of the particle in 
the direction perpendicular to the basal face). 




Fig. 4: Impact positions of random walkers from our simu- 
lations onto the surface of a hexagonal column, as viewed in 
projection from the side. The flux of walkers is highest on 
the sharp edges and corners of the particle, where the vapor 
density gradient is largest. 

closely approximates the data points (to within 1%) 
and this curve is overlaid in figure 5 for comparison. 

It is common in modelling studies to approximate 
hexagonal prism type particles as spheroids or circu- 
lar cylinders. The capacitance for a spheroid with 
major and minor axes matched to the length and 
maximal basal span of the hexagonal crystal is shown 
in figure 5, and there is reasonable agreement over the 
range of aspect ratios considered here to within 15%. 
Similarly, Smythe (1962) calculated the capacitance 
of a circular cylinder with diameter 2acyi and length 
L as: 

C = 0.637(1 + 0.868.A°-^'')aeyi. (4) 

It is interesting to note the similar form of equations 

3 and 4. According to Smythe, (4) is accurate to 
within 0.2% for aspect ratios between ^ = and 
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Aspect ratio = L/2a 

Fig. 5: Capacitance of hexagonal columns and plates with as- 
pect ratios between A = 0.01 (thin plates) and yl = 50 (thin 
columns). Diamonds are data points from our calculations in 
units of a. Solid line is fitted curve C/a = 0.58(1 + 0.95yi"'^"'). 
Dashed line shows the fit suggested by Chiruta and Wang 
(2005). Dotted curve is the capacitance of a spheroid with 
major and minor axes chosen to match the length L and basal 
span 2a of the crystal. Grey line is the capacitance of a circu- 
lar cylinder of radius a and length L as calculated by Smyth 
(1962). Cross indicates the metal-model result from the exper- 
iment of Podzimek (1966). 

^ = 8. Comparison shows that the hexagonal prism 
data Hes between the capacitance of a circumscribed 
cyhndcr (diameter 2a, length V) and an inscribed 
cylinder (diameter v^o, length L), as one would ex- 
pect since these cases consitute rigid upper and lower 
bounds for C (see appendix A). The curve for the cir- 
cumscribed cylinder is plotted in figure 5. 

Also plotted on figure 5 is the fitted curve derived 
for hexagonal columns and plates by Chiruta and 
Wang (2005) using a finite-element method. Their 
linear fit (7 = (0.751 + 0.49iyi)a is 10 25% higher 
than our data points over the range of aspect ra- 
tios which they considered (yi = 0.2-3.33), indicat- 
ing that their results arc overestimates. Their fit is 
also higher than that for an enclosing circular cylin- 
der, confirming that their data points are erroneously 
high (sec appendix A). We note that although their 
capacitances arc overestimates, Chiruta and Wang's 
main conclusion that solid and hollow columns have 
almost identical c;apacitances to one another is still 
likely to be correct, and this is in keeping with metal 
model measurements (Podzimek 1966). 

The metal model experiments of McDonald (1963) 
and Podzimek (1966) allowed estimates of the capac- 
itance to be made for hexagonal columns and plates. 
However there were significant sources of error in 



these experiments, and the results should be treated 
with care. Podzimek used an electrolytic tank to es- 
timate the capacitance of a metal hexagonal column 
with L = 50mm and 2a = 19mm (corresponding to 
A = 2.63), and measured that C = 19.5mm, ie. a 
normalised capacitance of C/a = 2.05 in apparent 
agreement with Chiruta and Wang's results. Unfor- 
tunately, Podzimek reported a systematic bias of ap- 
proximately -1-20% in his experiments (this was esti- 
mated by comparing the measured value for a thin 
circular disc with the theoretical value from table 1). 
To resolve this problem, he measured the capacitance 
for a metal model of a spheroid with the same width 
and length, and calculated the ratio of the measured 
column and spheroid capacitances. This ratio (1.116) 
was then assumed to be the same as the ratio of the 
true capacitances, and using the formula in table 1 
he deduced that C/a = 1.67. This point is plotted 
on figure 5, and is in good agreement with our fitted 
curve (C/a = 1.72). 

Similar experimental difficulties were encountered 
by McDonald (1963) who placed his metal models 
inside a walk-in Faraday cage and estimated the ca- 
pacitance of the arrangement. The connecting lead 
from the model to the capacitance meter shorted out 
many of the field lines which should have led from 
the Faraday cage to the model. This led to C being 
underestimated by as much as 45% for a thin circu- 
lar disc. His resolution of this problem was the same 
as that of Podzimek. and the ratio of the measured 
capacitance of the non-smooth models was estimated 
relative to idealised shapes. Prom these measured ra- 
tios, McDonald estimated that a hexagonal plate has 
a similar capacitance to that of an equivalent-area cir- 
cular disc. This is consistent with the data presented 
here: our fit (3) predicts the capacitance of a hexag- 
onal plate of zero thickness to be C = 0.58a, whilst 
an equal-area circular disc has the almost identical 
capacitance of C = 0.579a. 

6. 'Scalene ' columns and plates 

It has been observed (Bailey and Hallett 2004) 

that the basal faces of columns and plates are not 
always perfect regular hexagons, but are often some- 
what distorted. Here the capacitance of two so-called 
'scalene' forms similar to those reported in that paper 
is calculated, and the impact of breaking the perfect 
hexagonal symmetry is assessed. The shapes of the 
model scalene basal faces considered here are shown 
in figure 6. In both cases the maximum span across 
the basal face is defined as 2a, and the aspect ratio A 
is defined as before. The capacitances of these scalene 
crystals have been calculated for three aspect ratios; 
.A = 0.1, = 1 and A = 10, and the results com- 
pared with a perfect hexagonal crystal of the same 
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Fig. 6: Geometry of basal faces for the distorted or 'scalene' 
hexagonal columns and plates. For both types the maximum 
span across the basal face is defined as 2a and all of the inter- 
nal angles are 120°. The model shapes are: (a) a 'fiattened' 
hexagon with 4 short sides and 2 longer sides; (b) a second 
scalene type with 3 short sides and 3 longer sides. 

maximum basal span and length. 

The basal face of the first scalene crystal is a 'flat- 
tened' hexagon, with 4 short sides and 2 longer sides, 
and a width perpendicular to the maximum basal 
span of a (compared to \/3a for a regular hexagon). 
The internal angles arc all 120°. For A = 0.1, C is 
reduced by approximately 15% relative to a regular 
hexagonal plate; for ^ = 1 the reduction is 10%, 
and for a column with yi = 10 it is only 5%. 

The second model scalene crystal type has 3 short 
sides and 3 long sides, with a span between opposite 
sides of \/3a (the same as for a regular hexagon) . For 
this type it was found that C is essentially identical to 
a regular column/plate for all three aspect ratios, to 
within 3%. We conclude that distortion of the hexag- 
onal geometry of the crystals has a relatively small 
impact on their capacitance, and the only significant 
difference in the growth/evaporation rates is likely to 
be from surface migration and molecular accomoda- 
tion effects which are not included here. 

c. Bullets and bullet-rosettes 

Bullets and bullet-rosettes are one of the most 
common crystal habits in cirrus clouds (Heymsfield 
and laquinta 2000). Here we consider some simple 
model bullets and bullet-rosettes as illustrated in fig- 
ure 7. The bullets are modelled by a hexagonal col- 
umn of length L, width 2a, with a hexagon-based 
pyramid attached to one end. The ratio of the height 
of the pyramid cap relative to the column length is 
denoted by P, and in what follows we will assume a 
value of P = i. The aspect ratio of the columnar 
section of the bullet A = L/2a, is defined as before. 

The capacitance of single bullets with different 




Fig. 7: Model geometry for: a) single bullet, b) 2-rosette, c) 
4-rosette and d) 6-rosette crystal types. 

thicknesses corresponding to aspect ratios A between 
1 and 10 was calculated. The measured capacitances 
are plotted in figure 8. Equation 3 for hexagonal 
columns of the same width (= 2a) and total length 
{= L + PL) is also plotted on the figure. Comparison 
between the two shows that our model bullets have a 
capacitance which is reduced by a uniform 10% com- 
pared to that of a column. Reducing the size of the 
pyramid cap to P ~ ^ (following Macke et al 1996) it 
was found the capacitance is only 5% lower than that 
the complete column, and further reductions in P re- 
sult in an asymptotic approach to the capacitance of 
a simple column. 

Three model rosettes were constructed, and these 
are also illustrated in figure 7. The '2-rosette' is sim- 
ply a linear combination of the two bullets, joined at 
the tip. Calculations for a range of aspect ratios in 
the range = 1 to 10 indicate that the capacitance 
of the 2-rosette is reduced by 10-15% relative to a 
solid column of the same overall length and breadth, 
which seems reasonable given its geometry. 

The '4-rosette' has 4 bullets lying in a plane, each 
neighbouring bullet separated by 90°. Due to the 
more complex geometry of this crystal shape, its ca- 
pacitance has been normalized relative to the maxi- 
mum dimension D^ax (rather than the bullet width). 
Calculations have been performed for a variety of as- 
pect ratios: these values are shown in figure 9. The 
curve: 

C — 0.35yl " Dmaxi (5) 

was fitted to the data, and appears to provide a close 
approximation to it. Also shown in figure 9 is the ca- 
pacitance of an oblate spheroid with major and mi- 
nor axes matched to the overall dimensions of the 
rosette (major axis = 2L + 2LP, minor axis = 2a) 
and this curve overestimates the capacitance quite 
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Fig. 8: Capacitance of single bullet crystals as a function of 
aspect ratio. Diamonds are data points from our calculations 
in units of a. Dashed line shows equation 3 for hexagonal 
columns of the same width and total length as the bullet; solid 
line is that curve reduced by 10%. 

severely except when the aspect ratio of the bullets is 
close to unity. The capacitance for a circular disc and 
a sphere of equal maximum dimension are also shown 
on the figure for comparison. The value for a disc is 
a close approximation for A ~ 1.5, but increasingly 
overestimates C as ^1 increases. The capacitance of 
the sphere is a strong overestimate, typically a factor 
of two larger than the rosette. 

The final model crystal is a '6-rosette' made of 
six bullets in a three-dimensional cross shape, each 
seperated from its neighbours by an angle of 90° . The 
data for this rosette is also shown in figure 9 and the 
capacitances are approximately 15% higher than for 
the 4-bullet rosette. A curve was fitted to this data, 
and it was found that: 

C = 0.40.A-°-25i?„,, (6) 

is a good approximation to it. Again, a sphere of 
the same maximum dimension strongly overestimates 
the capacitance at all aspect ratios; a disc provides 
a closer approximation, but still overestimates the 
capacitance at large aspect ratios (thin arms). 

The results above for four- and six-arm rosettes 
with yi ~ 1.5 are consistent with the idealised rosettes 
modelled by Chiruta and Wang (2003). They used a 
series of smooth lobes to represent each bullet and 
showed that rosettes with four and six arms had a 
capacitance broadly similar to that of a circular disc 
with the same maximum dimension. An important 
feature of the new results presented here is that the 



Fig. 9: Capacitance of bullet-rosette crystals as a function of 
the aspect ratio of the arms A. = L/2a. The capacitances are 
normalised relative to the maximum dimension of the rosette 
Dmax- Diamonds are the sampled values of C for the 4-arm 
rosette; crosses are the values for the 6-arm rosette. Solid black 
and gray lines are the respective fitted curves (see text). The 
dashed line indicates C for a sphere of the same maximum 
dimension; the dotted line shows value for a thin disc. The 
capacitance of an oblate spheroid with a major axis of 2L+2LP 
and minor axis of 2a is shown by the dash-dot line. 

capacitance of bullet rosettes is sensitive not only 
to the number of arms, but also to the width of 
those arms. The aspect ratio of the bullets making 
up real ice rosettes (as observed by in-situ imaging, 
eg. Heymsfield et al 2002) is often rather larger than 
A = \.b and this means that the value of C may in 
fact be somewhat lower than that of a disc (~ 20% 
less for a four-arm rosette with A = 'i). 

d. Stellar and Dendrite crystals 

Two model stellar/dendrite crystal types were 
constructed. The first was a simple six-armed star 
shape, where each arm is a 'flattened' hexagonal plate 
similar to those described in section (3b). Each of 
these arms is separated from its neighbours by an 
angle of 60°. The arms overlap at the centre of the 
crystal, and the internal angles of the arms are all 
120°. The span between the tips of two opposite 
arms is defined as 2a, and the thickness is fixed as 
(2a)/100. The width of the arms is characterised by 
the ratio A' = w/a which we define as the ratio of 
the separation between the two longest sides of the 
hexagon w (ie. the arm width) to the arm length 
a. These dimensions are marked on figure 10, which 
shows the distribution of impacting walkers for a star 
crystal with A' = 0.2. 
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Fig. 10: Distribution of walkers impacting on the surface of 
star-type crystal (top) and dendrite with secondary branches 
(bottom). 

The capacitance of these star crystals is plotted as 
a function of the width of the arms in figure 11. For 
a star with thin arms {A' = 0.02) the capacitance is 
40% lower than a hexagonal plate of the same overall 
dimensions (= 0.596a for a plate of thickness 2a/100). 
As the arms become thicker, the capacitance rises and 
asymptotically approaches the value for a solid plate. 
The curve: 

C = 0.596(1 - 0.386-" ''^' )a, (7) 

is also plotted in figure 11 and approximates the data 

to within a few percent. 

The second model is an adaptation of the star 
model above, where a pair of secondary branches 
has been added to each arm, each with one third 
the dimensions of the main arms. These secondary 
branches are positioned with one end at the centre 
of the main arm, and oriented at an angle of 60° on 
either side of it. Figure 10 also shows the distribu- 
tion of impacting walkers on this crystal type, with 
A' = 0.05. 

The capacitance of this model crystal was calcu- 
lated for different arm thicknesses, and the data is 
shown alongside the results for the star shapes in fig- 



Ql , \ , \ , 1 , 1 

0.1 0.2 0.3 0.4 

Arm aspect ratio A' 

Fig. 11: Capacitance of stcUar/dcndritc crystal types. Di- 
amonds are for simple star shapes, crosses are for branched 
types. The capacitance of a solid hexagonal plate of the same 
overall dimensions is indicated by the dotted line. Dashed and 
solid lines are fitted curves for star and branched crystal data 
respectively (see text). 

ure 11. Adding the branches increases the capac- 
itance somewhat for dendrites with thin arms (for 
A' = 0.05 the branched dendrite has a capacitance 
~ 15% larger than a star). As the branches become 
thicker they approach the value for a solid plate. The 
curve: 

C = 0.596(1 - 0.30c-" **-^')a, (8) 

is a reasonable fit to the data, and this function ap- 
proaches the solid plate limit somewhat faster than 
the simple star crystal. 

e. Aggregates 

In this section we consider aggregates of the above 
ice crystal types. Aggregates are often the domi- 
nant particle habit in deep non-precipitating cirrus 
clouds (Field and Heymsfield 2003, West brook et al 
2006), as well as in snowstorms (Jiusto and Weick- 
mann 1973). A recent theoretical model of ice crys- 
tal aggregation (Westbrook et al 2004) has allowed 
us to produce large samples of realistic 'synthetic' 
ice aggregates. We have calculated the capacitance 
of these synthetic snowflakes in the expectation that 
the results ought to be a good approximation to the 
capacitance of natural ice aggregates. 

The random walker method was applied to ~ 1000 
synthetic aggregates sampled from the simulations of 
Westbrook et al (2004): a few examples are shown 
in figure 12. We have calculated the average capac- 
itance of the synthetic aggregates, using 10^ walkers 




Fig. 12: Examples of aggregates made up of 2,4,10 and 32 
hexagonal columns, sampled from the simulations of West- 
brook et al (2004). The aspect ratio of the columns is A = 4. 



to sample each individual aggregate. We note that 
this is a smaller number of walkers than was used 
for the pristine ice crystal habits; however, the com- 
putational cost is increased for complex shapes, and 
we arc averaging the results over many realisations of 
the aggregate geometry (every aggregate snowflake is 
different), so we expect that our eventual statistics 
should be accurate. In any case the results for the 
unit cube in scc;tion 2 indicate that the error in the 
calculated values for each individual aggregate should 
be less than 10%. 

For each aggregate, the ratio of capacitance to 
maximum dimension C/D^„^ was calculated, and the 
results binned and averaged as a function of the num- 
ber of crystals in the aggregate: this is shown in figure 
13. The 'monomer' crystals in this case were hexag- 
onal columns with an aspect ratio of ^ = 2. The ca- 
pacitance of the monomer columns is C / Dmax = 0.34 
(see section 3a), but for an aggregate of just two 
columns this is reduced to 0.21 ± 0.02, since two 
cohimns stuck together is a mtich more open geome- 
try than a single column. As more columns are aggre- 
gated the normalised capacitance rises to an asymp- 
totic limit of 0.25 ± 0.02. This value is in strong 
agreement with recent in-situ aircraft measurements 
of aggregate snowflake sublimation made by Field et 
al (2007) during a Lagrangian descent through a sub- 
saturated portion of an ice cloud, where a value of 
C/Dj^ax = 0.26 was estimated. 

The calculations described above were repeated 
for different monomer aspect ratios, and this is shown 
in figure 14. We find that for all the monomer shapes 
considered the normalised capacitance approaches 




# monomers per aggregate 



Fig. 13: Capacitance of aggregate snowflakes as sampled from 
the simulations of Wcstbrook et al (2004). The monomer crys- 
tals were hexagonal columns with A. = 2. The calculated val- 
ues of C/Dmax were binned and averaged as a function of the 
number of crystals per aggregate, as shown on the horizon- 
tal axis. Error bars are one standard deviation, points with 
no error bar indicate the value for a single aggregate. The 
normalised capacitance asymptotically approaches a value of 
C/Dmax — 0.25. Overlaid are theoretical values for a sphere 
(dashed line) and a thin disc (dotted line) with the same max- 
imum dimension. 

an asymptotic value of between 0.25 (thin columns 
A = 8) and 0.28 (squat columns A = 1). The area 
ratio of these aggregates (the ratio of the particle's 
projected area to the area of a circle of diameter 
Dmax) was measured and was also found to approach 
an asymptotic value of between 0.1 {A = 8) and 0.35 

(yi = i). ^ 

Changing the monomer crystal type to bullet- 
rosettes rather than columns is found to have very 
little effect on the capacitance, again yielding asymp- 
totic values in the same range. This capacitance 
is smaller than a disc or sphere of the same maxi- 
mum dimension, and the implication is that numer- 
ical models which assume these simple shapes are 
overestimating the growth/evaporation rate, by a fac- 
tor of two in the case of the sphere. 

It may perhaps appear counter-intuitive that the 
capacitance of aggregates should approach a con- 
stant value relative to their maximum dimension 
{C/Dmax = constant). The structure of the aggre- 
gates approaches a fractal geometry (Westbrook et 
al 2004) and as a result becomes increasingly open 
as the aggregates grow to contain more and more 
ice crystals. Because of this, the effective density 
decreases with size as {Dmax)~^ (in agreement with 
experimental data: Brown and Francis 1995, Heyms- 
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Fig. 14: The capacitance of aggregates with different monomer 
ratios A between 1 and 8. All of the curves approach an asymp- 
totic value in the range 0.25-0.28. 



field (it al 2002). One might therefore anticipate that 
C/Dmax would be reduced for large aggregates, since 
there is more empty space within a radius of 
the particle centre, and therefore (one would imag- 
ine) more opportunity for a water molecule coming 
within that radius to escape. However, the results 
of Ball and Witten (1984) show that water molecules 
following a Brownian path are exceptionally efficient 
at exploring three dimensional space, and because of 
this the ice particle appears essentially opaque to the 
incident water molecules (i.e. a fixed fraction of those 
venturing within D^axf^ will be absorbed). As a re- 
sult, C/D^ax = constant is in fact the physically sen- 
sible result for aggregates (see appendix B for further 
details). 



4. Screening 

Growth rates derived from laboratory experi- 
ments such as those of Bailey and Hallett (2004) 
are usually measured from crystals which grow not 
in isolation, but surrounded by other growing crys- 
tals. Screening of one ice crystal by another may be 
an important effect in light of the above discussion 
on aggregates, since the diffusing water molecules are 
very efficient at exploring the space around the parti- 
cle, so neighbouring crystals may constrict the vapor 
supply to one another. An understanding of such 
screening effects is important if laboratory data are 
to be accurately interpreted. 

For our model setup we consider 8 identical hexag- 
onal columns with ^ = 5 in a 'spiral staircase' geom- 



etry, growing with their longest axis in the horizontal 
direction. Neighboring crystals are offset in the ver- 
tical direction by a column's width and are rotated 
around the vertical axis by 45°. This is shown in fig- 
ure 15 and is roughly modelled after photographs of 
the experimental setup in Bailey and Hallett (2004). 
The surface vapor density and temperature are as- 
sumed to be the same on each crystal, and 10" walkers 
were used to sample the capacitance of the crystals. 
It was found that the vapor supply to the innermost 
crystals is inhibited to such an extent that their ca- 
pacitance is reduced to approximately one third of 
its value for a single column in isolation (C/a ~ 0.75 
compared to 2.41 in isolation). This result demon- 
strates how sensitive the vapor field around a crys- 
tal is to other sources/sinks of vapor in the vicinity. 
It also shows that if the growth rates from experi- 
ments are to be compared to theoretical capacitances, 
the crysals must be grown in isolation, or separated 
from one another as much as possible. Electrody- 
namic trapping techniques (eg. Swanson et al 1999) 
where ice particles are levitated in an electric field 
and allowed to grow in isolation may be useful in this 
respect, although to the authors' knowledge the re- 
sults have so far been limited to small particles, less 
than lOO/xm in size. 

5. Discussion 

The capacitance for pristine ice crystals and ag- 
gregate snowflakes has been calculated using a new 
Monte Carlo method, and these new values should be 
an improvement on the traditional approximations of 
smooth spheroids or discs. Application of these ca- 
pacitances to estimate the actual growth/evaporation 
rates requires a prescription for the difference be- 
tween the surface and the far-field vapor densities 
{ps — poo)- Using the Clausius-Clapeyron equation, 
an expression for dm /At in terms of the supersatura- 
tion with respect to ice {S — 1) is obtained: 



dm ^ ^ 5-1 

— — = A-kC X 

dt A + B 



(9) 



where the parameters A{T) and i3(T, P) are given 
in Pruppacher and Klett (1997) and depend on the 
ambient temperature T and pressure P. We believe 
that the results described in this paper may help to 
improve the estimation of growth and evaporation 
rates calculated in numerical weather prediction and 
cloud resolving models, especially for spatially ex- 
tended particles such as bullet rosettes, dendrites and 
aggregates, where modelling the particle as a simple 
sphere of the same maximum dimension (eg. Liu et 
al 2003, Khain and Sednev 1996) can overestimate 
|dm/di| by a factor of 2. This work may be particu- 
larly valuable to new precipitation models which pre- 
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Fig. 15: Model 'laboratory' setup of 8 identical hexagonal 
column crystals. Each is rotated by an angle of 45° about the 
vertical axis relative to its neighbours, and offset in the vertical 
direction by a crystal's width, in an attempt to model crystals 
growing on a vertical filament (indicated by the dotted line). 
The growth rate of the innermost crystals is approximately one 
third the growth rate of a column growing in isolation. 

diet ice particle habit from the model temperature 
and humidity (eg. Woods et al 2006), allowing an 
appropriate capacitance to be applied for each pre- 
dicted habit. We note that models which assume the 
capacitance of an equivalent volume sphere can lead 
to errors in C which arc size dependent, depending 
on the mass-dimension relationship assumed. Using 
equivalent volume spheres is an unphysical approach 
to estimating the capacitance, since C is determined 
by the physical dimensions of the particle and not 
its volume^. Parameterising C in terms of D^^x is a 
more natural approach since it represents the overall 
dimension of the particle and is also the parameter 
usually estimated from aircraft observations. 

The key limitation of the present study is the 
question of whether the growth and evaporation rates 
of ice crystals are simply dominated by the rate at 
which vapor impinges on/diffuses away from the ice 
crystal surface, or whether surface effects plays a 
significant role in limiting the growth/evaporation 
rates. For evaporation of ice particles, Nelson (1998) 
showed that the mass loss rate is controlled by (9) 



■^We note that the capacitance of an equal volume sphere 
is always an underestimate for the capacitance of any non- 
spherical shape (Polya and Szego 1951). 



and that surface migration of water molecules can 
be neglected, although molecular accomodation may 
play a role for very small, cold crystals (Magce et al 
2006). This conclusion appears to be supported by 
in-situ and laboratory imaging of sublimating ice par- 
ticles (Korolev and Isaac 2004, Swanson 1999) which 
indicate that the vapor flux away from the particles 
is concentrated at the corners and edges. This in 
keeping with the expectation for evaporation through 
bulk diffusion of vapor as described by equation 1 
without any surface migration. We therefore expect 
that equation 9 should be directly applicable to the 
problem of calculating the evaporation rates of most 
natural ice particles using the estimates for C pre- 
sented here. 

For deposition the situation is less clear cut: the 
appearance of the wide variety of different crystal 
habits at different temperatures and supersaturation 
levels indicates the influence of molecular accomoda- 
tion and surface migration on the growth. Despite 
this, the evidence from a number of expcrinicnits is 
that equation 9 is a reasonable approximation for the 
growth rate in many cases, although there is still 
much uncertainty as to the influence of molecular 
accomodation on growth rate at different tempera- 
tures and crystal sizes as summarised in Magee et 
al (2006), Fukuta and Takahashi (1999) and Prup- 
pacher and Klett (1997); in particular, poor accomo- 
dation of impinging vapour may be rather important 
for small, cold crystals. For frozen drops, laboratory 
measurements (Korolev et al 2004) have shown good 
agreement with (9). The new capacitance results pre- 
sented in this paper should allow a more accurate 
comparison to be made between theory and experi- 
mental data than has previously been possible, and 
may help experimentalists to estimate the 'accomo- 
dation coefficient' of complex ice crystals. 

Bailey and Hallett (2004) compared the growth 
rates of laboratory-grown plate and column crystals 
grown over a temperature range of —20 to — 70°C 
to the growth rate calculated using (9), modelling 
their lab crystals as spheroids (which our results indi- 
cate should be a reasonable approximation to within 
15%). Their results indicated that the predicted 
growth rates were only consistent with the measured 
values for aspect ratios close to yi ~ 1: for more 
extreme aspect ratios the theoretical values were as 
much as a factor of four too large for thin columns, 
and a factor of eight too large for thin plates. This 
may imply that the water molecules cannot be eas- 
ily accommodated at these low temperatures. Bailey 
and Hallett argue that because of the discrepancy 
between predicted and measured growth rates, the 
whole 'electrostatic' approach is unsuitable for cal- 
culating deposition rates at cold temperatures, and 
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they recommend that laboratory-measured growth 
rates should be used instead. However, the sensitiv- 
ity of (1) to screening as discussed in section 4 show 
that their measured growth rates may be strongly af- 
fected by neighbouring crystals, restricting the vapor 
supply, and reducing the measured capacitance sig- 
nificantly (a factor of 3 in our model setup). This 
may explain, at least in part, why their estimated 
growth rates are so much lower than equation 9 pre- 
dicts. The supporting glass filament may also play 
a screening role, particularly for thin columns where 
the filament thickness (50-70/um) is the same size or 
larger than the basal faces of the crystals. On the 
other hand, the thermal conductivity of this filament 
is much higher than the surrounding air, enhancing 
heat conduction and increasing the growth rate, lead- 
ing to further uncertainties. The interpretation of 
laboratory measured growth rates, and their use in 
evaluating the accuracy of theoretical models there- 
fore requires a great deal of care, and this is an issue 
where more experiments and theoretical work are ur- 
gently required to inform such comparisons. 

A key limitation of capacitance theory is that 
it is unable to make any prediction about particle 
habit, and a particular shape must be assumed a 
priori. Growth under equation 1 with the (moving) 
boundary condition p = constant at the surface is a 
much studied problem in theoretical physics: if the 
water molecules are simply deposited according to 
the distribution of vapor flux over the particle sur- 
face (eg. figure 4) the growth is unstable (Langer 
1980, Mullins and Sekerka 1963) and fern-like frac- 
tal patterns emerge (eg. Witten and Sander 1981, 
Bowler and Ball 2005; for a general review see Sander 
2000). That ice crystals grow in a much more con- 
trolled way indicates the infiuence that surface dif- 
fusion and the anisotropic molecular accommodation 
of the impinging water vapor have on the particle 
growth. It seems likely that dendrite crystals are a 
case in point, with their fern-like structure indicative 
of unstable growth, and the broad six-fold symmetry 
indicating the influence of the underlying crystalline 
anisotropy. This has recently been modelled using 
a crude anisotropy parameter by Goold et al (2005) 
to produce fractal dendrites with hexagonal symme- 
try. Progress in producing theoretical models which 
predict habit may well follow this approach, using 
quantitative data for the anisotropic molecular acco- 
modation coefficient and surface migration effects at 
different temperatures and supersaturations. 

The capacitance of realistic ice aggregates has 
been calculated for the first time (to the authors' 
knowledge), and the results are in close agreement 
with those estimated by Field et al (2007) from in- 
situ aircraft observations of evaporating ice aggre- 



gates, again indicating that (9) is a good approxi- 
mation for sublimating ice particles. It should be 
noted that the correction for the enhancement of the 
evaporation due to ventilation was estimated using 
formulas for idealised shapes (Hall and Pruppacher 
1976) and this introduces some uncertainty into the 
comparison. However, since a more accurate estima- 
tion of the ventilation effect would require a detailed 
treatment of the flow pattern around the aggregates, 
this value is at present the most accurate estimate of 
C that it is currently possible to make from aircraft 
observations. 

The authors believe that the capacitances pre- 
sented here are, in any case, the best estimates cur- 
rently available for the growth and evaporation of 
realistic ice crystals and snowflakes within the frame- 
work of the electrostatic analogy. 



Appendix A 

Here we show that the capacitance of a hexago- 
nal column with width 2a and length L (see figure 
3) must be less than that of a circular cylinder of 
diameter 2a and length L. 

Let the surface vapor density be = on the 
surface of the hexagonal column. Consider the flux 
of vapor through a cylindrical surface S'cyi with the 
above dimensions, just enclosing the column. The 
net flux of vapor $net through this surface must be 
the same as that onto the surface of the column since 
there are no other sinks of vapor present within the 
enclosed volume. 

As figure 16 illustrates, the flux can be split into 
two components: $n<,t = $cyi - ^escaped- The flux 
of water molecules incident from the outer boundary 
onto S'cyi for the first time is simply given by the flux 
onto a perfectly absorbing cylinder of the same di- 
mensions $cyi = 47r£)p<3oCcyi (since this corresponds 
to every water molecule trajectory being terminated 
at the point it first intersects S'cyi). The value of 
Ccyi is given by equation 4. Subtracted from this is 
the flux of water molecules which pass through Scyi 
but are not absorbed by the hexagonal prism and es- 
cape to infinity = ^escaped- Since this latter fiux is 
finite, $net < ^cyi, and by Gauss's law C < C^yi, 
ie. the capacitance of a hexagonal prism is lower 
than that of an enclosing circular cylinder. By the 
same argument it follows that a circular cylinder of 
diameter ^/3a must have a lower capacitance than 
the hexagonal prism which it inscribes. This allows 
the construction of upper and lower bounds on the 
capacitance of a hexagonal prism using the results of 
Smythe (1962). The violation of the upper bound by 
Chiruta and Wang's (2005) results shows that their 
data are overestimates. 
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Fig. 16: Cylindrical surface (dashed line) just enclosing a 
hexagonal prism (solid line), viewed in projection from one 
end. A trajectory which passes through the cylindrical sur- 
face and is absorbed on the hexagonal prism is illustrated; also 
shown is a trajectory which passes through the cylindrical sur- 
fcLce but escapes again and is not absorbed. 

Appendix B 

Here we show that C/D^ax is a function only of 
the shape of an ice aggregate, and not its size. Imag- 
ine that the particle is removed and it is replaced 
with a 'ghost' particle which is completely transpar- 
ent to the diffusing water molecules. The average 
water molecule follows a Brownian trajectory, trac- 
ing a fractal path with dimension d^, = 2 (Falconer 
2003). This means that whilst the water molecule is 
within a distance Dmax/'^ of the ghost particle centre, 
the resulting path fills a volume oc {Dj^axY™ ■ The 
ghost particle on the other hand, occupies a volume 
oc {DraaxY^^ With rf^ ~ 2 for our aggregates. 

Since + di > 3, the points where the particle 
and the uninterrupted random walk intersect with 
one another occupy a volume cx {Dmax)'^- Because of 
this, the aggregates appear essentially opaque to the 
diffusing water molecules, in the sense that for a given 
monomer type, a fixed fraction of the molecules which 
venture within a radius Dmax/'2 will come into con- 
tact with the ice particle. The ratio C/D^ax simply 
represents this fraction, and is therefore independent 
of size. 
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